#!/usr/bin/python
# -*- coding: utf-8 -*-

# Copyright (c) 2011
#
# Permission is hereby granted, free of charge, to any person obtaining a
# copy of this software and associated documentation files (the "Software"),
# to deal in the Software without restriction, including without limitation
# the rights to use, copy, modify, merge, publish, distribute, sublicense,
# and/or sell copies of the Software, and to permit persons to whom the
# Software is furnished to do so, subject to the following conditions:
#
# The above copyright notice and this permission notice shall be included in
# all copies or substantial portions of the Software.
#
# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
# SOFTWARE.
#
# Author: Jesus Carrero <j.o.carrero@gmail.com>
# Mountain View, CA
#

from scipy import linspace, zeros, arange
from scipy.sparse.construct import spdiags
from scipy.sparse.linalg import spsolve as spsolve

from TimeStepper import TimeStepper

import time
import logging

class zeroSecondDerivativeTimeStepper(TimeStepper):

    __solts__ = ['m_matStruct', 'm_MFwd', 'm_MBwd',
                'm_dataM', 'm_dataS', 'm_diags']
    """
        m_matStruct ~ 3, 5 diag or something else.
        m_MFwd ~ matrix to fwd step
        m_MBwd ~ matrix to bwd step
    """
    def __init__(self, method='cn'):
        TimeStepper.__init__(self, method)

        (self.m_MFwd, self.m_MBwd) = (None, None)
        (self.m_Ut, self.m_matStruct) = (None, None)
        (self.m_dataM, self.m_dataS, self.m_diags) = (None, None, None)

    def set_mat_attrib(self, struct):
        self.m_matStruct = struct

    def get_mat_structure(self):
        return self.m_matStruct

    def set_moperators(self, dataM, dataS, diags):
        (self.m_dataM, self.m_dataS, self.m_diags) = (dataM, dataS,
                diags)
    def set_boperators(self, MB, KB):
        (self.m_MB, self.m_KB) = (MB, KB)

    def setInitialSolution(self, u0):
        self.m_u0 = u0

    def get_history(self):
        return self.m_Ut

    def _solve_fwd(self, rhs):
        return spsolve(self.m_MFwd, rhs)

    def __theta_stepper(self, Utmp, MB, KB, i):
        """ Fwd step, update Blv and Brv using the rule:
            V_0^{k+1} + V_0^k = 2 * V_1^k  || V_m^{k+1} + V_m^k = 2 * V_{m-1}^k
        """

        Utmp = self.m_MBwd * Utmp
        (Blv, Brv) = (self.m_Blv, self.m_Brv)
        Blv[i] = 2 * Utmp[0] - Blv[i-1]
        Blr[i] = 2 * Utmp[-1] - Blr[i-1]

        theta = self.get_theta()

        Utmp[0] += (Blv[i - 1] - Blv[i]) * MB[0, 0] + (theta * Blv[i
                - 1] + (1 - theta) * Blv[i]) * KB[0, 0]
        Utmp[-1] += (Brv[i - 1] - Brv[i]) * MB[0, 1] + (theta * Brv[i
                - 1] + (1 - theta) * Brv[i]) * KB[0, 1]

        if 2 == MB[:, 0].size:
            Utmp[1] += (Blv[i - 1] - Blv[i]) * MB[1, 0] + (theta
                    * Blv[i - 1] + (1 - theta) * Blv[i]) * KB[1, 0]
            Utmp[-2] += (Brv[i - 1] - Brv[i]) * MB[1, 1] + (theta
                    * Brv[i - 1] + (1 - theta) * Brv[i]) * KB[1, 1]

        return self._solve_fwd(Utmp)

    def _init_engine(self):
        pass

    def step(self):

        # bwd, fwd matrices.




